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ABSTRACT 


An investigation was conducted to determine the 
feasibility of obtaining the short period and phugcid 
stabii ity derivatives from one maneuver, 
FouMuLrcaveousiLy"., it was concluded that the maxinun 
likelihood identification program SCIDNT I _ showed 
great promise in obtaining the short period and 
phugoid Sit a Dia iy: derivatives from one 
maneuver. Extraction of the short period stability 
parameters in the presence of the phugoid was 
easy, straight forward, and yielded results Similar to 
those obtained from pure short period data. Estimation 
of the phugoid stability parameters was possible when 
they were estimated in conjunction with 
Z, Xk , M , and ®. It was reconmended that anew set 
of sate ae obtained at a flight condition where the 
phugoid is at least moderately damped and that this 


data be anaiyzed to resolve the present anomalies. 
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A. BACKGROUND 


Flight stesting of military aircraft has always been a 
time consuming and expensive proposition. The relatively 
recent introduction of the new Paying Qualities 
Specification (#IL-F~&7&5B) has introduced many new 
parametric reguirements which require that many of the 
aircraft airframe parameters be determined to insure 
specification compliance. As a result of these new 
requirements the Naval Air Test Center (NATC) initiated a 
research progran to develop airframe parameter 
identification technology for use in flight testing Navy 
aircraft. This program waS initiated in 19771 between the 
Naval Airc Systems Command (NAVAIR) and NATC. At the same 
time the Office of Naval Research (ONR) issued a contract to 
Systems Control, Inc. (SCI) to advance the state~of-the-art 
in parameter identification. in 1973 a concentrated effort 
was begun by the Navy to continue the development of this 
new technology by forming a joint program between 
NATC, NAVAIR, ONR and SCI. AS a result of this development a 
digital FORTRAN conputer program, SCIDNT-I (hereafter 
refered to as SCIDNT) was written for NATC to extract 
aircraft stability and control coefficients from flight test 
data, Ref. 1. It was hoped that by using this computer 
program it would be possible to significantly reduce the 


fiight test time in the air, and the ground data reduction. 
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Be PURPOSE 


The purpose of this research was to develop parameter 
jdentification test methods for determining aircraft phugoid 
and short period characteristics from one aircraft 
maneuver. This research further entailed the determination 
of proper aircraft maneuvers/pilot inputs and detailed data 
analysis proceedures using the advanced maxinum likelihood 


identification technigue. 


C. DESCRIPTION OF ALRPLANE 


The test airplane was demonstration F-14A ship number 8 
BU. NO. 157967. The airplane was essentially typical of a 
production F-14A airplane relative to these tests. An 
aircraft description is contained in Appendix A. Photographs 
and a three view drawing of the test airplane are shown in 


Appendix A , Figures 1, 2 and 3. 


D. GAXIUUN LIKELIHOOD PARAMETER IDENTIFICATION THEORY 


A widely used rule for choosing the sample statistic to 
be used in parameter estimation is called the method of 
maximum likelihood Ref. 2. This method, in additicn to 
prescribing the statistic which shouid be used, provides an 
approximation of the distribution of that statistic so that 
approximate confidence intervals may be constructed. The 
nethod of maximum likelihood makes use of the 
sample-likelihood function. Assume we observe or record a 
Set of observations on a system, whose model has p unknown 
parameters 6. We then can write the joint probability 


distribution of a random sample for any given set of values 
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of the parameters @ from a feasible set ®. Re then can 
assign a probability p(zj{6) to each outcome 2z. If we 
observe or measure the outcome of an actual experiment to be 
z we then would like to know which set or sets of values of 
@ might have resulted in these observations. This concept is 
contained in the likelihood function © (@[z). This function 
is of fundamental importance in estimation theory because of 
the likelihood principle of Fisher and others which states 
that if the system model is correct, all information about 
unknown parameters is contained in the likelihood 
function. The maximum likelihood estimator of ®@ is the value 
@ which causes the likelihood function & (O®{[z) to be a 


maximum [2]. 


6 = max ££ (6{2z) 
Ge @ 


From a purely mechanical viewpoint this method can be 


summarized as follows: 


Pind the probability density functions of the 
observations or measured vaiues for all possible 


combinations of the unknown parameter 
vaiues. Select the density function whose value is 
hignest among ali apa | fines lLous at the 
measured values of he observations. The 


corresponding parameter values are the maxinunm 
likelihood eStinmates. 
Assune ® can take three possible values: O oer and 
: Further assume that the probability density functions of 


observations z for these three values of 6 are as shown in 


Figure 1. Then, if the actual observation is Z, S is the 


maximum likelihood estimate of 9. 
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observed 


Figure 1. Maximum Likelihood Estimates 


Owing to the monotonic, one-to-cne relationship between 
the likelihood runction and its logarithm, i.e. they have a 
maximum at the same value 8, it is more convient te work 
With the negative of the logarithm of the likelihood 
function. Starting from a priori values, the parameters are 
updated so that the value of the negative log-likelihood 
function J at the observed value of the outputs decreases 
continuously. 

As Sample statistics le functions of the 
Observations, the estimators are, before the experiment, 
random variables and can be studied as such. In fact, the 
major advantage of maximum likelihood estimators is that 
their properties have been thoroughly studied and are well 
Known, virtually independently of the particular density 
function under consideration. Additionally, the maximum 
likelihood method can be used with linear or nonlinear 
models in the presence of process and meaSurement noise. The 


estimator of 6 1S approximately normally distributed for 
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large n, and is asymptotically unbiased, efficient, 
consistent, and sufficient, [2]. 
The details of the maximum likelihood identification 


proceedure are presented in the following section. 


E. MAXIMUM LIKELIHOOD IDENTIFICATION OF PARAMETERS 


The derivation follows closely that presented in Ref. 
3, and 1S presented herein for completeness and reader 
convenience. The initial equations will be concerned with 
general ncnlinear systems. | | 


Consider the general nonlinear aircraft equations of 


motion 
mee it (x, U, @, t) + | (e, t)w Oma =< 7 (1) 
E(x(0)) = x (8) 
O 
i 
oii) - x (67) (2 (0) - x (07) } = P (8) ( 2 ) 
O O O 


where 


Ce) 1s n x 1 state vector 

u(t) iSgeeex 1 anput vector 

6 1S p x 1 vector of unknown parameters 

["(8, t) is n X q process noise distribution ntatrix 
w(t) 1S q x 1 random process noise vector 


Sets of m measurements Ye) are taken at discrete times t 
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y(t) = BSCS Me BSG 6, on + ee (3) 


k = 1, a Dace. Sieg N 


w(t) and ac) are Gaussian random noises [2], with the 


following precperties 


E(w(t)) = 0 B(v(t,)) =a) 


Tr 
E (w(t) v Ceo 0 


Jb 
E(w(t) w (7)) Q(6,t) 6 (t-T) 


T 
evga) ¥ (taj) = R(S,t ) co, oo 
3) k 5 ak 


The unknown farameters are assumed to occur in the 


functions f and h and in matrices [". Gee he Po. and “X.« 2h 
O O 


the following analysis, the model and the functional form of 


f and h is assumed known correctly. 


The set of observations Nahe (lo), «seo py (tC) 
1 2 N 
constitutes the outcome Z in this case. The likelihood 


function for 6, which has the same form as the probability 
of the cutccme z for a certain value of parameters 6, is 
given by 


x (8/2) = p(z{9) 


PAVE De V(t )reeee¥(t ) 18) 
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L(Ofz) = Saat (5) 


where 


Y= (y(t deren eV (t de K = Ve 2 eee 


Y {@) = ee) 9,0) © ply 9 
p(Y 16) = p(y(t iy eye PW 18) 


, 9 sf e 
) P(X 18) 


= c Y , © t Et 
Piy(t )IY 0, ©) piy(t, 4¥ 


which by successive application of Bayes rule, [2], becomes 
N 
p(Y {8} = | p(y(t.)'%. + 9) ( 6 ) 
N vee a 2 ea 
1L1= 


The log-likelihood function is 


N 


log (XH(61z)) = ». log MEMO AGE IS 32 6)} + constant (7 ) 
a 


@O find the probability distribution of y(t_) given Y. and 
7 2 


8, tne mean value and covariance are determined first. 


E(y(t.)tY. ., 0) QO ¥(ifi-t) (8) 
pe 1-1 


The expected value or the mean is the best possible estimate 


of measurements at a point given the measurements up to and 
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including the previous point. 


“A a S A ° ° a 
orev t(c JY 976) = Ei{y(t.) = y(iJi-1)) (y(t) - y(iti-1)) } 
1 1-7 1 SE 


de 
QE{V (i) Wi) } 


Ly B (i) (293) 


v (i) are called the innovations at point i and B(i) is the 


innovations covariance. Since 


y(t.) - E(y(t.)1Y¥. .,0) = V (i) ( 10 ) 
a 1 i-1 


it follows that y(t.) given Y. and @® have the sane 
aL ie 


distribution as WEN It has been shown that as the 


Sampling rate is increased, the innovations V (i) tend 
towards having a GausSian density. AsSSuming a sufficiently 
high Sampling rate, the distribution of V (i) and, 


Eherefore, y(t.) given Y. ; and © is Gaussian, i.e., 
1 : 


1 T 
exp|- 3) B-1 (i) va | 
Peyqe ) 1Y. ,98) = 
ri i-] 


/ {2 


m/2 west 
(2 Tr) |B (2) | 


( 71) 


we 





1 T 
log {p(y(t.)1¥. _-6)} = -—V (i) B-1(i) Vi) 
¥ 1-1 Z 


1 
> Logie (s) | + constant 


( 12 ) 
The log~likelihood function of equation ( 7 ) can be written 
as 
N 
1 | an 
log (£(81z)) = mre [v~ (4) B12 (4) V (4) 
+ logl!B (i) 1} (oa 


An eStimate of the unknown parameters 1S obtained by 
Maximizing the likelihood function or the log-likelihood 


function from the reasible set of parameter values. 


®> 
tl 


aan log ( £(81{2Z)) ( 14 ) 
1 : T 
= nax | - =) {fy (i) B-1(i) V(i) + log{B(i) 0 | 


2 a 
6Ea@ 1=1 


( 15 ) 


The log-likelihood function depends on the innovations 
and their covariance. To optimize the likelihood function, a 
way must ke found for determining these quantities. Both 
innovations and their covariance are outputs of an extended 


Kalman filter. 
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The extended Kalman filter is conventionally divided 
into two parts. In the first part, cailed the prediction 
equations, the state equations and state estimate covariance 
equations are propagated in time from one meaSurement point 
to the next. In the Second part, called the measurement 
update equations, the measurements and associated 
Measurement noise covariances are used to improve state 
estimates. The covariance matrix 1S also updated at this 
point to reflect the additional information obtained fron 


the measurements. 
1. Prediction Equations 


The state prediction is done using the equations of 


motion, Appendix 8B. Starting at time t. With current 
2 ge 


estimate %(i-1[i-1) of the state x(t.) and the covariance 
ae 
P(i-iji-1), the following ecuations are used to find the 


predicted state xXx(iji-1) and the associated covariance 
P(iji-1). 


d aA a - 
ames ie = r] Ve 
a x ee) A ee UGE) 8, t) ( lo ) 


. T 
Poeepeoe ) - CAC) ePYttt ) et P(tit ) F {(t) 
1-1 . i-7 1-1 


ay 
+ Poof ( 17) 
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The nx n matrix F is obtained by linearizing f 


about the best current estimate 


Ae (xu dt) ,u(t),6,t) 
1-1 


ee — 18 
a gx (tit. _) 
1-1 


Using ( 16 )-( 18 ), we can obtain 


x(t it, .) Qe (iti-1) 
2 14-1 


and 


P(t {t ) QP(ili-1) ( 19 ) 
eel tect 


Thereafter, the measurement update eguations are used. 


2. Measurement Update Equations 


ee ee eee eee 


The covariance and state estimate are updated using 
| the measurements. The necessary relations are derived in 
Ref. 4, and are presented here without proof. The innovation 


and itS covariance are 


v(i) = y(i) - h(X(i1i-1),u(t.),@,t,) ( 20 ) 
a Te ae 


| | | ; 
B(i) = H(i) P(iJi-1) H (i) + R ( 21) 
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where H iS obtained by linearizing h, 


gh (x (i 4i-1) miley o,c ) 
1 L 


1 NS 22 
sa dx (i}i-1) 


The Kalman gain and the state update equations are 


T 
Key =P ii-?) A (i) Bq (i) ( 23) 
X(ii{i) = X(ifi-1) + K(i) VW (i) ( 24 ) 
femmad ly, © (14), the covariance of error in updated 


State, 1s obtained by 
P(ili) = (I-K(i) H(i)) P(iji-1) ( 25m) 


One can now return to the tine update or prediction 


eguations. 


we OPtrapizetion Procedure 


— ww se es Se ae a ee ae a ee Se ee ee oe 


Many possible numerical procedures can be used for 
this optimization problem. Modified Newton-Raphson or 
Quasilinearization have been found by experience to give 
quicker convergence than most procedures such as the 
conjugate gradient or the Davison method. The modified 
Newton-Raphson method is a second order gradient procedure 
requiring computation of first and second order partials of 


the log-likelihood function. 
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The last three terms in the equation for the second 
Pastial SimeeiceLog-likelihood function involve second 
partials of innovation and its covariance. These terms are 


usually dropped which yields the second partial 
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: approximation 
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From eguations 20 and 21 we can obtain the gradients 


of innovation and its covarience for parameter ® aS, 





a v (4) gh dx(iji-1) «gh 
ae x je ”:—is cin 
c 5 x=X (i{i-1) 
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dB (1) gH (1) ots i es (lee |) oT 
LI ~1)H H 4 
10 36 Pie it (1) + (2) 10. (2) 
J 5 J 
H i R 
+ H(i) > ( 30 ) 
J J 


j = eeigietals 2 
fy OS as eo a 


Recursive equations can be obtained for gradients of 


the predicted state and its 
stages by using the 
equations of the Kalman 


foeto with respect to @ yields 


covariance. 
prediction and measurement 


filter. Differentiating 


This is done in 


update 


equations 


J 
gad )CCgE (2 (t1t, ),0(t),t,8) 
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dt ge. ; de | 
J J 
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The sensitivity functions, used in the SCIDNT integration 
routines, are updated at measurement points by 


feELerentidemig equations 22 -~ 25 with respect to 6 . 
J 


2h(X(ili-1),u(t ),0,t 
gH (i) 7h (x (ifi-1) ,u( 2 al 


SO 34 
19 © g8 4x (it i-1) 
i, J 
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J = 1,2,3,-+-,7P 


The negative of the matrix of second partials of tae 
log-likelihood function is called the information matrix 


M. The step size A6© for parameter estimates is given by 


1 ) 
em 205 1 e812) ( 38 ) 
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The information matrix provides a lower bound on parameter 


estimate covariances, 1.@e., 


A ae 
E (6-6) (9-6) 2 M7! ( 39 ) 


This is the Cramer-Rao lower bound. The maximum likelihocd 


estimates approach this bound asymptotically. 
4. Linear Systems 


In a linear systen, the functions f£ and h are 


defined as 
f (x,u,9,t) QF(@,t)x + G(e,t)u 
and h (x,u,@,t) QH(68,t)x + D(8,t)u ( 40 ) 


For the linear case, it has been shown that 1£ the 
model is correct, the innovations are white and have 
Gaussian density at the true values of the parameters. The 
assumpticn of fast sampliag rate 1s, therefore, not 
necessary, 

The basic algorithm for solution of the linear case, 
is the same, however, some of the equations can now be 
Simplified. The equivalence between F(t) and H(i) of ( 18 ) 


and (22) #=.and F(@,t) and H(@,t) of (40) is 
ak 


obvious. Equations ( 16) and ({ 20 ) now become 
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o A 
— = - t ’ 
AL oS are F(9,t) x pee + G(®@,t) u(t) 

( 47) 


and 


Vv (i) = y(i) - H(O,t.) X(ili-1) - D(®,t.) u(t.) 
1 L i 


( 42 ) 
Equations ( 29 ) and ( 31) can be written as 
; Nee: qH(6,t_) 
d ¥ (i) dx (i{i-1) ‘fe 
UUM A X (71 1 1 1) 
d9 | 1 ¢9 | do | 
J J J 
dD(8,t.) 
ee (t (43) 
el 
4 
and 
dx(tit. .) gx(ti(t. ) 
i | th 
a (8st) 
dt dS | dO | 
a J 
F(9,t) , 
jy ail La ) 
(9 | a=" 4 
J 
G(8,t 
Uae u(t) ( 44 ) 
qo | 
J 


All otner equations remain the same. 

There iS a considerable reduction in the computation 
reguirement for time-invariant ilinear systems. In this 
case, Matrices F, G, H, D, [~, Q and R and their derivatives 


with respect to the parameters are constant. 
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MemivertanteLincar Systems In Statistical Steady 


a> —— — SO —— —— a ee 


In many aircraft applications, the Kalman filter is 
in steady state for the duration of the experiment. This 
occurs when the Kalman filter is in operation for a 
sufficiently long time and the process and measurement noise 
covariances are constant. The Kalman gain and the 
innovations and the state covariances approach constant 
values. The time update and measurement update equations for 


the covariances are 


: PCE i el (ae [ene Bele ; Q : 
sP(tit. |) (tit, (tit of 


( 45 ) 
Te 
K = P(iji-1)H B-? ( 46) 
T 
B = HP(iji-1)H + R ( 47 ) 
P(iji) = (I-KH) P(ifi-1) (48 ) 


By definition of the steady state 


P(i-1{i-1) = P(ijfi) 


Therefore, from ( 45 ) 
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QQ g( at) (I-KH) P(ifi-1) g(at) + Q! 
T 
= g(4t) (P(ili-1) - KBK ) g( at) + Q! 


CS) sy 


Using ( 46), ( 47 ) and ({ 49 ), we can solve for P(iji-1) 


and then find K and B. Also, it can be shown that 


oF A dP a i PA P 50 ) 
— = —_— + 4 = 
@ 286=6h1 qo. ae a 
I-KH aE rs HPH R)-1 + A 51) 
—_—_—_— = a —— a = oF 
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Thus, it is possible to solve for 





Wognmg( of ) and then 


dK gB : 
vs 7+ ao Erom  { 51 ) and ( 52). Equation ( 50 } 





is a linear eguation in and the coefficient of the 


unknown matrix does not opepend on the parameter 
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Serio, the Sensitivity of the state covariance matrix can 
J 


be determined very quickly for all parameters. Once the 
sensitivity of P, K and B for unknown parameters is 
determined, only state sensitivity equations need to be 
updated. The computation of state sensitivity functions can 
be reduced to many fewer equations. 

An approximation simplifies the problem further. The 
unknown parameters are defined to include elements in K and 
B matrices instead of Q and R. Optimizing the log-likelihood 


function for parameters in B gives 


fl 
eons 


N 
T 
Ye Yo (i) ( 54) 


The gradient of the log-likelihood function with respect to 


Other unknown parameters 1S 


N 

d Log( £(61[Z)) rT. a d V (i) 

aoa Se SS Vo(i) Bri —————— 55 

6 J (2) 0 ( ) 
j 1=1 y 


The sensitivity of innovations to parameters is determined 


uSing the following recursive equations 


da A : 
— x(t{t = Fx (t{t Gut 
eit) (tit, .) (t) 
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( 57 ) 
K(iji) = x(ifi-1) + K VW(ji) 
: F Ay io D 
dvi) sd ee ee. Oe 
{9 | qo | d& 9 | i 
2 2 J 3 
( 58 ) 
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gk . 
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Bhs, 





where I, , is a matrix of all zeroes except a 1 at the 
fe, kt Heetilion. | 
This approximation simplifies the optimization 
considerably. However, this usually leads to an 
overparameterized model. In other words, once K and 8B are 
determined, it is not possible to find any corresponding 
[". Q and R which have the desired (Known a _ priori) 
structure. A good estimate of elements in r, Q and R 
matrices can be obtained by a jeast-sgquares type 
approach. The fit to the observed data is better than with 
true values of parameters but the parameter estimates do not 
have minimum variance. This approximation is not good in 
aircraft application where the structure of Q and R is known 
fairly well, but is excellent where there are many process 
naise sources and the characteristics of both Q and R are 


relatively unknown (e.g., economic systems). 


ikelihood With No rocess Or Measurement 


The maximum likelihood method can be simplified when 


either precess noise Or measurement noise are 
absent. Assuming no process noise we can proceed as 
folloys. 


If the process noise is zero and initial states are 
known peptectly, i.e., W(t) and P (0) are zero, the 
covariance of the error in the predicted state is aliso 
zero. It is clear from ( 23) that the Kaiman gains are 


zero. The innovations are the output error, i.-e., 


Weljeasy (ajo h(x (t)),u(t ),0,t ) ( 60 ) 
ae at 1 
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and the innovation covariance 1s ( 21 ) 


B(i) = R ( 61) 


the log-likelihood function is, 


N 
1 fh 
log (3€(81z)) = a . Yo (i) BoA (i) + log 1Rf 
i-1 
(Gon) 


which on optimizing for unknown parameters in R yields 


N 
; 1 . Ta 
R ce Vo (i) ( 63 ) 


1 
ci 


The equality in ( 63 ) holds only for those elements of R 
Which are not known a priori. For instance, even if R& is 
known to be diagonal, the right hand side matrix will not be 
diagonal in general, but the off-diagonal terms should be 


ignored before they are equated to R. Using (63) in { 62 ) 
we obtain 


N 
y| fk 
log (#(@jfz)) = ->—) Vv (ye R= 2 NG} 
i=1 
+ constant ( 64 ) 


uy 





The optimizing function is the same as that for the output 
error method except that the measurement noise covariance 
matrix is determined using (63) and 1S used aS the 
weighting matrix in the criterion function. In the output 
error method, the measurement noise iS assumed known and the 
weighting function is arbitrary. 

The first and second derivatives of the 
log-likelihood function with respect to unknown parameters 


ale 


d i , ov (i) 
se 109K (812)) > yea) EE C65) 
7 1=17 5 
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The terms in the second derivative are approximated as 


: N ih 
g@log( £(01z)) y dv G) s , dv) 
aS) 8 8 S) 
co do a coe d j 
( 67 ) 
With no measurement nolse we can proceed as 


follows. If all states are measured with no noise, tne 


Covariance of the error in state estimates is zero at the 
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beginning Of any time update, 


tt 
or) 


P (i-1}i-1) 


( 68 ) 


X (1-1) 


and %(i-1{i-1) 


It is easy to show in this case that for fast sampling the 
log-likelihood function is quadratic in the difference 
between measure values of x and fie, tio tne method 
reduces to the equation error method, the weight W being 


chosen as 


a 
1 Ak 
W [= fo etc t) ) (x = £ (x,u,O,t)) dt 
ve) 
( e ) 
Thus, the maximum likelihood method and eguation-error 


methods are equivalent except for the technigue for choosing 


the weighting matrix. 


F. AIRCRAFT PARAMETER IDENTIFICATION USING THE MAXIMUM 
LIKELIHOOD METHOD 


In applying the maximum likelihood theory at is 
necessary tc consider the aircraft equations of motion in 
discrete form. It was assumed that the structure of the 
nodel was Known. The vector of unknown parameters from x(Q), 
Oo, G, RB poten, 9 Oe and R was denoted by 8. The 
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problem, therefore, was to identify 98 using the maximun 
likelihood method. In constructing the problem it was 
further assumed that the aircraft which had been flying for 
some time was ina steady state condition, and process noise 
(wind gust effects) was not present. 

In the computer implementation, the negative 
log-likelihood function was minimized using a Gauss-Newton 
gradient procedure. A least squares routine was utilized to 
provide approximate airfrane parameter Stant au 
values. Additionally the results of the SCIDNI computer 


program were utilized to obtain frequency and damping. 


fee SCOPE OF RESEARCH 


Research was limited tc the longitudinal case since both 
available computer and research tine were limited. Although 
additional data were desired, the project airplane was 
substantially damaged by fire and was removed from further 
test activity. Further data gathering with the bo ga 
airplane Was, therefore, not possible. Spee rac 
longitudinal stability derivatives identified are given in 
the Results and Discussion section. There were no airplane 
flight restrictions which affected the test methods or 
limited the scope of research. The research was limited to 
establishing base line data for the short period and phugoid 
modes and thence to determine the feasibility of identifying 


these modes Simultaneously. 
wee UETHOD OF TESTS 

In general test proceedures were in accordance with 
those contained in the U. S. Naval Test Pilot School Fixed 


Wing Stag it y and Control Manual, Ref 5. There 


were, nowever, some non-standard test techniques utilized in 
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exciting the desired airplane modes. The short period node 
was excited by the standard doublet method and additionally 
by applying a sinusoidal control input at varying 
frequencies. The phugoid mode was excited by the Delta 
Airspeed method in all cases, {5}. 

In exciting the phugoid and short period modes 
Simultaneoulsy, a combination of standard test technigues 
were utilized. The "simultaneous" technigue consisted of 
establishing the airplane in a phugoid maneuver and during 
the first quarter cycie of the maneuver applying a doublet 
input to excite the short period mode. 

The data system utilized within the airplane was 
specifically designed to be compatible with the Naval Air 
Test Center Real-Time Telemetry Processing System (RTPS). A 
brief description of the alrpiane and ground data systems is 
presented in Appendix C. Data reduction was performed cn the 
RTPS Xerox Sigma 9 computer in the batch processing 
mode. Raw airplane data wer2 converted to an Engineering 
Unit tape and used aS input to the digital parameter 
identification program (SCIDNT). In addition a least squares 
parameter identification routine was utilized to generate 
parameter start up values for SCIDNT. This least squares 
routine extracted a set of approximate stability derivatives 
by minimizing the error between the model and actual 


airplane state vector time histories. 
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A. SHORT PERIOD AIRFRAHE PARAMETER IDENTIFICATION 


Initial trimmed test conditions for each control 
input are presented in Table I. Prior to each control input 


the airplane was trimmed for steady state "hands off" flight 


[>]. 


Altitude MACH CG Loading 
(ft HPC) (% MAC) 
C1) 
Sine Wave Input Das ASIA. 26149 A118 B 
Doublet Input 30 503 nOecaeo O29 35 B 


(1). Four Sparrow III missiles (AIMN-7) on stations 3, 4, 5 
and 6 


Table I. Short Period Test Conditions 


The airplane model equations used were those shown 
in Appendix B. In order to isolate the short period mode 
from the longitudinal equations of motion the state vector 
parameters « , q, and 9 were nade insensitive to u. This was 


accomplished simply by fixing Xs Ae He 2 oe and Lee 
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equal to zero and not identifying any of them. This equation 


modification resulted in the clasical short period 
approximation Ref. 6. | 
TO provide SCIDNT with approximate airframe 
parameter start-up values,the least squares routine was 
utilized. For the short period case the least Squares Start 
up values were very good, a comparison of a typical set of 
least squareS parameters and the final SCIDNT parameter 


values is shown in Table II. 


2 vA M M Z M 


cc q a q Ce ores 
Least -.4909 0.96246 ~1.78 =e oa Ua 0 62/09 
Ssawares 
SCIDNT OS NOs Iau  —i.do <2 S/S. 02276 Bis ee 
Table Ii. Comparison of Least Squares Parameter 
Values and SCIDNT Parameter Values 
In actual practice, however, airplane data was 


recorded in long bursts which usually included three or more 
maneuvers , executed With near identical initial 
conditions. AS a result once a set of start up values were 
determined for a given test condition it was unnessary to 
determine a set of start up valves for other like 
maneuvers. 

The actual parameter identification proceedure , for 
the short period mode, was relatively Simple once 
familiarity with the digital program was attained. Upon 
obtaining a set of "reasonable" initial parameter 
estimates, all six short period parameters were estimated in 
an iterative fashion. Here "reasonable" denotes proper 
parameter sign and magnitude. Each computer run was made 
With the number of iterations and step cuts both initialized 


at ten. During a given iteration the program attempted to 
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minimize the likelihood function. If the likeiihood function 
increased, a series of step cuts would be taken until the 
likelihood function again decreased or the maximum number of 
step cuts were exceeded. The program would then continue to 
iterate until the desired number of iterations/stef cuts 
were reached. Due to program logic one nore step cut’ than 
requested was executed. For the short period 
mode, however, convergence was obtained early in  progran 
execution. In general, iterations conducted with doublet 
control input data converged very quickly and began 
executing step cuts until the allowable number was 
exceeded. Figure 2 shows the change in the likelihood 
function during iteration for a doublet input. Progran 
convergence had been essentially reached by the fourth 


iteration at which time a series of step cuts were made. 


Iteration 
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Figure 2. Likelihood Function Change During Iteration For 
Doublet Input 


The convergence characteristics using sinusoidal 


input data were somewhat different. During the progran 
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iteration the likelihood function continuously decreased and 
no step cuts were taken, Figure 3. Program convergence for 
this case was evidenced by a leveling out of the likelihood 
function, but not until the eighth iteration. It was felt 
that since the data length for the Sinusoidal input was 
longer than the short period data length the larger number 


of iterations until program convergence was not abnormal. 


Iteration 
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Figure 3. Likelihood Function Change During Iteration For 


Sinusoidal Input 


After convergence was obtained for the parameter 
estimates, the Angie Of Attack Scale Factor (K, ) was 
determined. Ki, was a scaling factor in the angle of attack 
measurement equation, Appendix B. in estimating i two 
terns in the model were actually affected, Se itself and 
the variable (1. Kx /U) which was used to multiply gq in the 
measurement equation. Since, in essence, two terms were 
being estimated a modified form of parameter eStimation was 
utilized. In the tern (ly K 70) the value of Ky was 


initially set to 1.0 and the single tern io was 
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determined. This new value of Ke was then substituted into 
the tern (ly Kec /U) and again a set of iterations were 
performed to determine an updated K,. This method of 
iteration was continued until there was essentially no 
change in Ki 7 at which time the value was said to have 
converged. In most cases only two iterations were required 
to obtain a converged Ky, . 


Upon determining a value for kK final set of 


7a 
iterations were - performed to Pe ebine values for 
Z,X ,M, 6. At this point, the short period mode time 
Pestory “fits were nearly identical and the estimate 
statistics indicated a high degree of confidence. Attempts 
at further iterations to improve the estimates resulted in 
using a lot more computer ese eb iess very little 
improvement. Therefore, the procedure utilized to estimate 
the short period mode was to first estimate the short period 
parameters, second estimate K , and finally, estimate the 


initial values. 


=e Se ae —_— 


Data from typical short period estimation runs for a 
doublet and sinusoidai input are presented in Tables III and 
IV. Typical short period time history fits are presented in 
Appendix D Figures 1 and 2, for the doublet and sinusoidal 
inputs respectfully. It should be noted that tne scales for 
the short period time histories are different for the 
doublet and Sinusoidal inputs. Furthermore, since the 
estimated stability derivatives were in dimensional form the 
actual derivative vaiues differed slightly from maneuver to 
maneuver. 

The data obtained with the doublet input in most all 
cases converged to believable derivative estinates, l.e., 
expected sign and magnitude, with favorable iteration 
statistics. For example, refering to Table III, the value of 
M was -1.844 with a standard deviation of 0.0061 and a 


x 
final step size of -0.0146 which is significantly lower in 
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magnitude than the estimated parameter. Fron probability 
theory we can say that the probability is 99.74% that the 
value of My is -1.844 +£+0.0183 which is well within 
acceptable engineering tolerances. This same analysis was 


carried out for the other estimated derivatives with similar 


results. 

Derivatives Z Z M M Los Ms 
a g eC g e e 
Estimate -.4138 1.099 -1.844 -+-.826 0.1628 5.098 
Standard 0.0024 0.0038 0.00651 0.0037 0.0055 0.0104 
Deviation 

Parameter ~.0006 ~.012 -.0146 -.0001 0.07143 9.028 
Step Size 

Eigenvalue ize ax 106 

Spreaa 

Likelihood 

Function -5.1 x 103 

Value 

Tine 

History O = 0.00396 OP = 0.564 @ = 90,0039 
Matches cS u q 


Gm Of 00s OC = 2.33 
‘s) n 


zZ 


Table III. SCIDNT Estimates of Longitudinal Short Period 
Stability Derivatives Using a Doublet Input 


The eigenvalue spread, that is, the difference 
between the largest eigenvalue and the smallest eigenvaiue 
of the information matrix, was 1.28 x 106 which was less 
than the maximum spread allowed by SCIDNT (1019). All 
parameters converged indicating no identifiability 
proplems, and the eigenvalue spread was therefore not 
considered significant. The likelihood function value was a 
large negative number which converged rapidly (Figure 


2), further indicating no identifiability problems. 
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As a final check of stability derivative correctness 
a series of time history fits were constructed. In 
constructing these plots both the model and airplane were 
excited with the same control input and the airplane tine 
history plots superimposed on the same scale. The tine 
history fit for the doublet input was, in general, very 
good, with the difference between the measured value and 
estimated value within the instrumentation measurement error 
tolerances Ref. 7. 

Referring to Table IV, it can be seen that typically 
the estimated derivatives using a sinusoidal input converged 


tc appropriate (sign and magnitude) stability derivative 


values. 

Derivatives Zi Z. M M Zo Nes 
a G c q S e 

Estimate =o te tea2o) =—2.007 =.7615 0.1532 5.2274 

Standard O20G42 02.0061 §0.013 0.007 0.0045 0.0181 

Deviation 


Parameter VO ea eOO To .004s =—20061tT —-.0037 0.0332 
Step Size 


Eigenvalue tGHerx 102 
Spread 


miKelihood 

Function -~7.323 x103 
Value 

Time 


History G = 0.00533 GF = 0.9 GT = 0.006 
Matches a u g 


Table IV. SCIDNT Estimates of Longitudinal Short Period 


Stability Deravitives Using A Sinusoidal Input 


in general, however, the parameter standard deviations and 
step sizes were not as favorable as those obtained with the 


doublet input for these two examples. It was 
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possible, however, to select another pair of maneuvers and 
reach the opposite conclusion, ie. Sinusoidal input more 
favorable. Since all parameter estimates for both the 
doublet and Sinusoidal inputs converged it was concluded 
that the slight difference in statistical data was not 
Peognieticant. Lf a determination as to which control input to 
use were required, it was felt that the primary 
consideration would be ease of pilot input. If an autopilot 
device were available which could input the sinusoid then 
this input should be used, however, with no autopilot systen 
installed the doublet input would provide excellent 
results. 

A consistent problem that occured throughout the 
Mescarch was concerned with the initial conditions prior to 
the control input. If the initial condition of the airplane 
was not stabilized prior to the control input or if a data 
burst was taken in the middle of a control input, no tine 
history match was possible. hecording to theory 1t is 
possible to extract the airframe parameters at any time 
during a maneuver, this, however did not appear to be 
true. Figure 3 Appendix D shows the results of eStimating 
the parameters starting from non-zero estate delete 
conditions. This time history fit was generated after ali 
iterations on both parameters and initial conditions were 
completed. The parameter estimates did not converge and the 
most powerful derivative, Moe was much larger in magnitude 
than any cther estimates with a value of -3.7. These results 
occured any time an attempt was made to start SCIDNT from an 
initial condition other than steady state. It was concluded 
penat, within the scope of this research, meaningful 
parameter estimates were only obtainable when started from 
an initial steady state condition. 

A comparison of Short period Characteristics 
obtained fron Navy Preliminary Evaluation (NPE) data and 


those determined by SCIDNT is presented in Table V. 
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Description Damping Frequency 


Rat1o 
C Ww) 
Sp n 
Sp 
(rad/sec) 
SCICEDNT-Doublet 0.403 1.54 
SCIDNT-Sinusoidal 0.371 1.61 


Table V. Comparison of Longitudinal 
Short Period Characteristics 


It was concluded that, within the scope of this 
research, the parameter identification program SCIDNT was 
easy and straightforward to use in determining the short 
period mode. In addition, it was concluded that it was not 
possible to obtain accurate estimates when starting from 


non-steady state initial conditions. 


Be. PHUGOID AIRPRAHE PARAMETER IDENTIFICATION 


i. Airplane Test Conditions 


a cial et ll agree 


Initial trinmed test conditions for the phugoid 
Maneuver are presented in Table VI. Prior to the centrol 
input the airplane was trimmed for steady state "hands off" 
meaght {5}. 


Altitude MACH CG Loading 
(Ft HPC) (% MAC) 
Phugoid Input 30 055 0.674 10.94 B 


Table VI. Phugoid Test Conditions 


34 








The airplane model equations used were those shown 


in Appendix B. The derivitives os and X» were fixed at 
2 
zero values since their effect on the phugoid was 


small. Additionally, it was felt, fixing these terms would 
facilitate the identification of the dominant phugoid 
terms. Finally from Appendix B it can be seen that the term 


Xx’ was defined as (X - @& U) and therefore due to its 
q -g 0 0 


relatively large size vaS approximated by cc JU and 
| oe) 


fixed. The parameters which were identified during the pure 
phugoid maneuver were therfore, Oe awe and nee 

Since there was no short period input during the 
pure phugoid maneuver it was not possible to effectively 
identify the short period parameters. Additionally the 
progran solution diverged without a set of predefined short 
period parameters and it was necessary to obtain parameters 
fron a prior run. A set of short period parameters 
were, therefore, obtained from a previous identification run 
and fixed throughout the phugoid parameter identification 
runs. 

In an attempt to provide SCIDNT with approximate 
airframe parameter start-up values the least squares routine 
was utilized. For the phugoid case the least squares start 
up values were unuSable. For example, least squares valves 


Or X was -2.6697 and Xa Was -~-240.5, both of which are 
e 
terms of small absolute value. As an additional check, the 


least squares parameter estimates were used to generate a 
set of time histories. For these time historys the standard 
deviation of airspeed and normal acceleration was 49.12 
ft/sec and 5.44 ft/sec® respectively. Therefore, since the 


least squares routine did not provide acceptable start-up 
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values the identified parameters were set initally to zero 
for the first SCIDNT iteration. Further attempts to obtain 
mogoid Start up values for SCIDNT from the least squares 
routine proved futile and its use was discontinued for the 
phugoid case. 

The version of SCIDNT utilized hou this 
investigation would accept a maximum of 500 data samples per 
run. Since the phugoid period wWwasS approximately 100 sec a 
data burst recorded at 20 samples per second (SPS) would 
only inciude the first 25 sec of the maneuver. This tine 
Span was obviously acceptable for the short period mode but 
clearly unacceptable during a phugoid maneuver. An attempt 
was therefore, made to utilize data at 2 and 5 SPS to 
include as much of the phugoid maneuver as possible. This 
data was from the same maneuver with data tapes recorded at 
meand 5 SPS. 

Although the 2 SPS data would provide a longer 
phugoid time span in theory, the program exhibited severe 
instabilities during the control input and therefore only 
the latter portion of the data burst was usable. It was 
strongly susvected that these instabilities were due to the 
unusually rceugh control inputs. 

The actual parameter identification procedure for 
the phugoid mode was somewnat different than that utilized 
during the short period investigation. The short period 
Stability derivatives were obtained from a run whose initial 
conditions most closely approximated those of the phugoid 
run. These short period parameters were then fixed as 
constants during the phugoid iterations. Computer runs were 
made with the number of iterations and step cuts both 
initialized at ten. The program would then continue to 
iterate until the desired number of iterations/sterp cuts 
were completed. For both the 5 SPS and 2 SPS data 
convergence waS obtained prior to executing the full 10 
iterations. Iterations conducted with phugoid data recorded 


at 2 SPS appeared t0O converge very quickly and began 
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executing step cuts until the allowable number was 
exceeded. AS in the short period case the program executed 
one step cut more than reguested. Figure 4. presents the 
change in the likelihood function during iteration using 
phugoid data recorded at 2 SPS. The program converged to its 
Mhest" solution by the third iteration at which time the 


step cuts were executed. 


Tteration 


-7700 





Likelihood Function Vaiue 


aod 9690999990000 


Pigure 4. Likelihood Function During Iteration 
Using 2 SPS 


Note that although the program had apparently 
converged at the third iteration the value of the likelihood 
function decreased until the fourth step cut. The function 
value then remained essentially constant throughout the 
remaining step cuts, 

The convergence characteristics using data recorded 
at 5 SPS were somewhat different, Figure 5. | 
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Figure 5. Likelihood Function During Iteration 
Using 5 SPS 
It should be noted, however, that the scales are 
different for each plot. The program had essentially 


converged at the fifth iteration, at which time the progran 
executed step cuts. Note that the Likelihood function 
values decrease uSing 5 SPS data, and were essentially the 
Same uSing the 2 SPS data. However, the time history matches 
for each set of data were significantly different. The angle 
of attack scale factor USS ) was not estimated since the 
value of Ke from the short period data was utilized. After 
convergence waS obtained for the parameter estimates a final 


set of iterations was performed to estimate values for 


Za Xt,’ , 8 and the three phugoid parameters 
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Simultaneously. At this point, the phugoid time history fits 


were acceptable for the 5 SPS data and clearly unacceptable 
using the 2 SPS data. Further attempts to improve the 
estimates and time histories resulted in using a lot more 


computer time with no improvement. 
3. Data Analvsis 


Data from phugoid estimation runs using 2 SPS and 5 
BES are presented in Tables Vaid and ones 
respectfully. Phugoid time history fits are presented in 
Appendix D Figures 4 and 5 for 2 and 5 SPS respectfully. It 
Should be noted that the scales of the phugoid time 
histories for 2 and 5 SPS are different. 

The data obtained using a Sample rate of 2 SPS was 
motally unacceptable. For example, refering to Table 


VII, the value of M was -0.00000017 while the last step 
u 


eaze was 0.002. 





Derivatives Z x X X M X 


u oc u q u Oe 
Estimate -.003 0 -~.158 wie) ~ ln 7xtO-?. 86 6¢ 
Standard -00006 fixed ~0542 ice -0002 fixed 
deviation 
Parameter -~.001 fixed 0.444 Fixed 0.003 fixed 
eee Pp Size 
Eigenvalue 3.54 x 108 
Spread 
Likelihood 
Function ~8.4 x 103 
Value 
Time 
History G0 = 0.0048 QO = 26.77 @ = 0.00534 
Matches oe u g 


Ga— 0.055. GC = 2.53 
.S] n 
Z 
Table VII. SCIDNT Estimates of Longitudinal Phugoid 


Stability Derivatives Using 2 SPS Data 


Although the eigenvalue spread and likelihood 
function change during iteration indicate no significant 
program difficulties the standard deviations of the state 
error vector are very large. AS an approximation it can be 
Said that the phugoid maneuver is essentially a varying 
airspeed and pitch attitude maneuver with a constant load 
factor. From Table ViI it can be seen that the standard 
deviation for airspeed was 26.77 f£t/sec-! and for pitch 
attitude was 0.053 rad, both of which were unacceptable. 

The time history plots constructed using the 2 SPS 
data snowed no real agreement between the model and the 
airplane. In constructing these plots both the model and 
airplane were excited with the same control input and the 


airplane time history plots superimposed on the Sane 


60 


» ¥ 


> 


oe] 





Beasve. £t 1S obvious that there was Simply no fit of any of 
the state vector. It should be noted that the plot of 
average horizontal stabilizer position although appearing to 
be very rough was a result of the plot scale only and was in 
essense relatively smooth. 

The parameter identification results using the 5 SPS 
data indicated it was possible to identify the phugoid mode 
parameters. Referring to Table VIII, it can be seen that the 


estimated derivatives using the 5 SPS data converged to 


appropriate (Sign and magnitude) Stability derivative 
values. 
Derivatives Z X X X M Xo 

u oc u g u e 
Estimate =, 00:0 2 0 ~.006 oso ono x tO > 20 


Standard eee V0m © > fixed -0002 mired wo. xX 10>? fixed 
Deviation 


Parameter -8.7x10-6 Fixed -.0016 fixed -.00002 fixed 
meep S1ze 


Eigenvalue 1.82 x 1032 
Ppt ead 


Likelihood 
Function -1.49 x 10% 
Value 


Tine 
History G@ = 0.00267 oO = 5.53 g = 0.00343 
Matches oc u g 


CGO = 0.014 gG~ = 0.92 
6 n 
Z 


Table VIII. SCIDNT Estimates of Longitudinal Phugoid 
Stability Derivatives using 5 SPS Data 


The standard deviations of the parameter estimates were much 


Smaller than the parameter, and the parameter step sizes 
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indicated that the estimates had converged to the "best" 
estimates obtainable with the data. It should be noted that 
standard deviations for u and 6 were 5.53 ft secy! and 0.014 
rad respectively, a significant reduction in standard 
deviation frem those obtained with 2 SPS data 
The time history plots using 5 SPS data _ showed 
relativeiy close agreement between the airplane and model 
time histories. There was, however, some model instability 
noted during the first 30 seconds of the time histories. It 
was concluded that these instabilities were due to the 
relatively slow sampling rate during the initial control 
mnput. 
From SCIDNT the phugoid damping ratio ( ¢) and 


frequency (w (rad secum!)) were 0.0486 and 0.06 
n 
respectfully. Available F-~14 data was obtained in the auto 


Sweep mode with SAS on while SCIDNT data was obtained in 
auto sweep, SAS off. It was therefore not possible to make 
any meaningful comparisons. However, it was noted that 
between 10,000 ft and 35,000 £t the phugoid period varied 
between 106 sec and 72 sec. The SCIDNT frequency resulted in 
a period of 92.53 sec which when considered in light of 
available NPE data appeared reasonable. 

It was concluded that, within the scope of this 
research, the paramater identification program ScipiNa 
provided a reasonable estimate of the phugoid longitudinal 


parameters. 


C. SIMULTANEOUS SHORT PERIOD AND PHUGOID PARAMETER 
IDENTIFICATION 


1. Airplane Test Conditions 


eae asa = ae SR a Pe See ae Se ee a a ee ee ce ee 


Initial trimmed test conditions for the short period 
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and phugoid maneuvers are presented in Table IX. Prior to 
the control input the airplane was trimmed for steady state 
"hands off" flight (5]. Since ths short period input was 


applied after the phugoid input, the two test conditions are 


different. 
Altitude MACH CG Loading 
(Ft HPC) (% MAC) 
Phugoid Input 30 326 ~6217 10.86 B 
Short Period Input 30 703. 2963 10.86 B 


Table IX. Short Period and Phugoid Test Conditions 


Zoe iciued OF Envesti gation 


Eee SS Se eS ae ES ae ee = 


The airplane model equations used were those shown 


in Appendix B. The derivatives xX , X,, and X were fixed 
ce On q 


as described in the previous section. 


The maneuver investigated in this section was 
designed specifically to provide the means of identifying 
the short period and phugoid parameters from one maneuver 
or, if you will, simultaneously. Initiation of the phugoid 
Maneuver was accomplished in accordance with standard test 
techniques, {5]. Upon establishing the airplane in the 
phugoid maneuver, a doublet input was applied to excite the 
Short period mode. The airplane was then allowed to fly 
through at least one phugoid cycle. 

Identification of the short period parameters was 
essentially routine and straightforward. The proceedures 
used in estimating the parameters were those descriked in 
the previous short period section. The computer runs were 
made with the number of iterations and step cuts both | 
initialized at ten. Although this short period input was 


Superimposed on the phugoid, convergence was obtained early 
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in program execution. Iterations conducted under these 
conditions converged very quickly and executed step cuts 
until the allowable number was exceeded. Likelihood function 
change during the short period iterations was similar to 
that shown in Figure 2. 

Based on the results of the previous Section it was 
concluded that data recorded at 5 SPS would most likely 
yield reasonable phugoid parameter estimates. The short 
period stability derivatives obtained from this run were 
fixed as constants during tne phugoid iterations. Computer 
runs were then made with the number of iterations and_ step 
cuts both initialized at ten. 

Figure 6 presents the likelihood function value 


during iteration of the phugoid stability derivatives. 


Iteration 
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Figure 6. Likelihood Function During Iteration 


Por Phugoid 
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The likelihood function change for the phugoid maneuver 
displayed some characteristics heretofore not 
encountered. The characteristics shown are typical of the 
results obtained with the available data. At some time 
during the iteration sequence the likelihood function would 
exhibit a jump discontinuity, in this case it occurred 
between the first and second iterations. Program iteration 
then appeared to proceed in a normal manner until the tenth 
iteration at which time the program terminated. From Figure 
6 it was obvious that convergence had not been obtained at 
program termination. 

In an attempt to correct the iteration instabilities 


the phugoid stability derivatives Z@, X , and WM were 
u u 


estimated Simultaneously with Z, X , M , and 6. Figure 7 
O fe) fe) re) 
presents the likelihood function during a typical iteration 


for seven parameters. 
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fieerati1on 


—3000 


Likelihood Function Value 





Jae job ok: J. Likelihood Function During Iteration For 


Phugoid (Stability Derivatives, 2, X , M, and 6) 
| Oo O Oo O 


It can be seen that the likelihood Function value converged 
rapidly to a steady state value with no apparent 


SEeSCOHEINULtLIeS or nonliaearaties. 


3. Data Analysis 


—_y ee oe —_ 


Data from a typical short period run superimposed on 
the phugoid are presented in Table X. A typical time history 


fit 1s presented in Appendix D, Figure 6. 
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Derivatives Z Z i M Zo Ms 
ao q oa g e e 

Estimate -.4285 1.048 -1.54 ~ .665 0.19 G 464 

Standard 0.0049 0.0074 0.014 0.011 0.0082 0.0312 


Deviation 


Parameter 0.0149 0.0119 -.0166 0.0355 0.0173 -.1815 
Step Size 


-Eigenvalue 1.25 x 105 

Spread 

Likelihood 

Function -~4.36 x 103 

Value 

Time 

History G7 = 0.00188 Gm O0nce72 FG = 0.0036 
Natches x u g 


G7 = 0.0048 Gg = 1.82 
8 n 
74 
Table X. SCIDNT Estimates of Longitudinal Short Period 


Stability Derivatives 


The short period parameter estimates displayed 
characteristics Similar to those observed when estimating 
the short period parameters from a short period input 
only. The time history plots showed good agreement between 
the model and actual airplane data. It was noted that 
airspeed was not constant since the short period was 
Superinposed on the phugoid. This varying airspeed did not 
however, have any noticable effect on the determination of 
the short period parameters. 

The results for the phugoid portion of the maneuver 
were at times in disagreement with theory. From theory 1t 
can be shown that the speed damping derivative, which acts 
along the negative X axisS increaseS with speed and is 
therefore negative [6]. During the iteration to determine 


the phugoid stability derivatives it was noted that the 
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speed damping derivative (X ), became positive on the 
u 
seventh iteration and remained positive until program 


termination, Table XI. 


Derivatives Z xX X X M X 


u ° u q u oO 


e 
Estimate -.0002 0 ~065 onc —~seoxto-S& 0 


Standard Verte10-© fixed 0002 fixed 3.5x10-© fixed 
Weviation 


Parameter 5.38x10-© fixed 4.8x10-5 fixed 1.2x10-5 fixed 
weep S1ze 


Eigenvalue 6.2 x 1012 
Spread 


Likelihood 
Function -1.005 x 10% 
Value 


Time 
History T= 0.004 GF = 1.97 G = 0.0047 
Matches a u gq 


G = 0.031 G = 1.48 
6 n 
Zz 
Table XI. SCIDNT Estimates of Phugoid Stability 


Derivatives (Positive X ) 
u 


The tine history plots were in reasonable 


agreement, Appendix D Figure 7. The positive value of X was 
u 


not however considered reasonable since this value would 


imply decreasing drag with increasing airspeed. The X tern 
u 


was however negative through the sixth program iteration. As 


mentioned previously the likelihood function during 


iteration showed definite signs of instability, it was 
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therefore suspected that the speed damping Sign change was 
due to thas program instability. The program was 
therefore, run again and terminated at the sixth iteration 


thereby producing a negative speed damping term, Table XII. 


Derivatives Z X X X M Xo 
u Oo u g u e 
Estimate - .00008 0 =~.Q00032 -43.85 0.00017 0 


Standard 2-6x10—-& fixed 0.00047 fixed 6.7x107S fixed 
Deviation 


Parameter 4.0x10-5S fixed -.0059 fixed 0.000058 fixed 
Strep Size 


Eigenvalue 8.53 x 1011 

Spread 

Likelihood 

Punction =f Cee 104 

Value 

Time 

History 7 = 0.0071 G-~ = 3.94 G =0.0053 
Matches oc u gq 


Oo = 0.0366 @ = 2.14 
6 n 


Zz 


Table XII. SCIDNT Estimates of Phugoid Stability 


Derivatives (Negative X ) 
u 


With the iterations thus restricted to six the value for X 
u 


WaS maintained at a negative value. Time history plots for 
the negative X case are presented in Appendix D, Figure 
u 


8. It was noted that there was no appreciable improvement in 


the time history flots. In fact from Tables XI and XII it 


can be seen that the standard deviation of the state error 
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vectors was degraded for the negative X case. This 
u 
condition in fact is however, not very startling when it is 


realized that the program had demonstrated instabilities and 
had not converged. 

In an effort to reconcile these inconsistencies and 
lack of convergence the three phugoid parameters were 


estimated along with 2, X , M, and @. This’ technique 
O O O O 


yielded reasonable stability derivatives and acceptable time 


MstOby ElOCS. Stabality deriyative data provided by the 
seven parameter identification are Shown in Tabie XIII, time 


history rflots are presented in Appendix D, Figure 9. 


Derivatives Z X x X M Xo 


u ee u q u & 


Estimate =e OO Ore O eeu  aNI6G> =. 00014 1 0 


Standard 1.4x10-& fixed MU comertrxede 5.2xX10=¢ fixed 
Deviation 


Parameter 1.23x10-% fixed 6.3x10-5 fixed ~-4.02x10-© fixed 
Beep Size 


Eigenvalue 3.7 x 1033 
Spread 


Likelihood 
Function ~4.23 x 103 
Value 


Time 
History T= 0.0045 QO = 4.17 gq =0.00523 
Natches x u g 


G =0.031 T = 1.86 
8 n 
Z 


Table XIII. SCIDNT Estimates of Phugoid Stability 


Derivatives (Parameters Estimated With Z , X , M , and ® ) 
O O ) O 
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The parameter standard deviations were smail and 
Within engineering tolerance. The parameter step sizes were 
quite small and showed further evidence of solution 
convergence. The standard deviation of the calculated model 
response with respect to the measured airplane response was 
acceptable. The tine history plots provided further 
indication of parameter estimate validity. 

The problems associated with the phugoid using the 
Simultaneous technique could have been a result of progran 
deficiencies or airplane anomalieS or a combination of 
both. The flight conditions at which the data was obtained 
were areas in which the phugoid mode was lightly damped or 
neutrally damped. In addition the airplane wing was set to 
the auto sweep mode. It was felt that either of these 
conditions couid have caused the program instabilities 
noted. Attempts to obtain new data were not possible due to 
airplane non-availability as mentioned earlier. 

It was concluded that extraction of the short period 
stability parameters in the presence of the phugoid was 
easy, straight forward, and yielded results similar to those 
obtained from pure short period data. Estimation of the 
phugoid stability parameters was possible when they were 


estimated in conjunction with Z,X , WM, andé. 
fo) O O re) 


It was recommended that a new set of data be 
obtained at a flight condition where the phugoid is at least 
moderately danped and that this data be analyzed to resolve 


the present anomalies. 


71 





Pte CONCLUSTONS 


a ae Eg eg 


A. GENERAL 


The maximum likelihood identification program SCIDNT 
showed great promise in obtaining the short period and 
phugoid stability derivatives from one maneuver, 


"Simultaneously". 


Be oO PeCITRFIC 


Wat n ip the scope of this research, the parameter 
identification program SCIDNT was easy and straight-fcrward 
to use in determining the short period mode. 

Meaningful parameter estimates were obtainable only when 
started from an initial steady state condition. 

Within the scope of this research, the parameter 
identification program SCIDNT provided a reasonable estimate 
of the phugoid longitudinal parameters. 

Extraction of the short period stability parameters in 
the presence of the phugoid was easy, straight-forward, and 
yielded results similar to those obtained from pure short 
period data. 

Estimation of the paugoid stability parameters was 
possible when they were estimated in conjunction with 


Z i X v My, and 8 ° 
oO O oO oO 
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IV. RECOMMENDATIONS 


= ee aes ee ame sm a ee 


It is recommended that a new set of data be obtained at 
a flight condition where the phugoid is at least moderately 
damped and that these data be analyzed to resolve the 


present anonalies. 
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APPENDIX A 


ERO PAR GRALT DESCRIPTION 


The F-14 airplane is a supersonic, two-place, 
twin-engine, swWing-wing, air-superiority fighter designed 
and manufactured by Grumman Aerospace Corporaticn. In 
addition to its primary fighter role, carrying missiles and 
an internal 20-millimeter gun, the airplane is designed for 
fleet air defense and ground attack missions. Photographs of 
the test airplane are shown in Figures 1 and 2. A three-view 
drawing is presented in Figure 3. 

The forward fuselage, containing the crew and electronic 
eguipment, projects forward from the mid fuselage and wing 
glove. Outboard pivots in the highly swept wing glove 
support the movable wing panels, which incorporate integral 
fuel cells and fuli-span leading edge slats and trailing 
edge flaps for supplemental lift control. In flight, the 
Wings may be varied in Sweep, area, camber, and aspect ratio 
by selection of any leading-edge sweep angle between 20° and 
689. Wing sweep can be automatically or manually controlled 
to optimize performance and thereby enhance aircraft 
versatility. Separate variable-geometry air inlets, offset 
from the fuselage in the glove, direct primary airflow to 
cwo waro0=P-U12A, dual-axiail-compressor, turbo-fan engines 
equipped with afterburners for thrust augmentation. The 
displaced engine nacelles extend rearward to the tail 
section, Supporting the twin vertical tails, horizontal 
tails, and ventral fins. The middle and aft fuselage, which 
contains the main fuel celis, tapers off in depth to the 
rear where it accommodates the speed brake surface and 


arresting hook. Retractable vanes in the glove leading edges 
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extend to supplement lift and compensate for changes in the 
Lr Grart aerodynamic center. All control surfaces are 
positioned by irreversible hydraulic actuators to provide 
desired control effectiveness throughout the flight 
envelope. Stability augmentation features in the flight 
control system enhance flight characteristics and provide a 
more stable and maneuverable weapons delivery platform. fhe 
tricycle-type, forward-retracting landing gear is designed 
for nose gear catapult launch and carrier landings. Missiles 
and external stores are carried from eight hardpoint 
stations on the center fuselage between the nacelles and 
under the nacelles and wing glove; no stores are Carried on 
the movable portion of the Wing. The fuel systen 
incorporates both inflight and Single-point ground refueling 
capabilities. 

A more complete description of the F-14A airplane is 


given in Ref. 8. 
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3/4 FRONT VIEW 


FRONT VIEW 


F-14A AITFPLANE 
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APPENDIX B 


DERIVATION OF THE AIRPLANE MODEL AND MEASUREMENT EQUATIONS 


A detailed derivation of the complete aircraft equations 
of motion is shown in Ref. 6. A Summary of the assumptions 
made during that derivation is shown below for 
completeness. These assumptions lead to a set of force and 


moment equations. 


* Summary of Derivation ASSumptions * 


1. The airframe is a rigid body 

2. The earth is fixed in space 

3. The earth's atmosphere is fixed with respect to the 
earth 

YW. The mass of the airplane remained consStant for the 
duration of any particular dynamic analysis 

5. The xz plane is a plane of symmetry 

6. The disturbances from the steady fiight condition is 
Small enough so that the products and sguares of the changes 
in velocities are negligible in comparison with the changes 
thenselves 

7. The disturbance angles are small enough sucn that the 
Sines of these angles may be set equal to the angles and the 
cosines set equal to one. Products of these angles are also 
approximately zero and can be neglected 

8. Since the disturbances were small, the change in air 
density encountered by the airplane during any disturbance 


can pe considered zero 
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9. During the steady flight condition, the airplane is 
flying with wings level and with all components of velocity 


zero except U and W 
O O 


10. Body axis system 1S assumed 


Dividing the force equations by the mass m and the 
moment equations by the appropriate moments of inertia 


yields terms of the forn: 


onl - and 1 ers (eB) 
u I qr 


1 

Ti 
Replacing the first term by X and the second term py L 
Simplifies the notation. These . quantities are called 
"dimensional stability derivatives". Using this short hand 
notation tke equations of motion applicable to the 
longitudinal problem can be written as foliows: (a complete 


set of eguations is presented in Ref.5) 


Ge 


+ Wgq+g@cosO® =Xut X.u+Xgqg*+txgq 
O O u a q G 


meeree wt ke So + ig C 4+ Xy O 
W W ae ae =: oo = 


+ X + Xn CO + Xn GO + XB CG 
Op O; Op “F O53 B Diy as 


ee S tcos Tu + cosé 2 
Cp Ws “uy Cann O ppy 1D j 
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Seeueget geusingo = Zu +t 2.8 + Z 4 + Z,q 
fe) oO a q 


ee FZ Ot 25 6 + 24 6 
W W eS eo © eye 


+22 6 + Ze 6 + Ze G 
F Op “Pp On “3 Tr 


+) 754 - } Tu-~ (si T a) B3 
ce (sin € ) u in § ! Onpm RPM 


51) 


get eet Meg +e g + Mmwet Mk + Ne OG 
u u gq gq W W e e 


+ Ms 6 +E GO + Me O + Mao O + Me OC 
ee once op ep Op OR OE 


Zz Mm Zz Mm 





ie ,° + Ng +———T u + ( Bu ) 


fh re) ne: 
B B L u i RPM RPM 
BOE yy 


In order to simplify these equations it was assumed that 
all disturbances were small such that the higher order terms 
were negligible. Also it was assumed that the flaps and 
speedbrakes were retracted and that thrust effects were 


negligible. With these further assumptions the w equation 


became: 
w - O q + g6 sin@ =Zu+t+Zqgt aw + L6 oO ( BS ) 
O ‘e) u q W e e 
let was further assuned that the flow was 


quasisteady. Because of this all derivatives with respect to 
the rates of changes of velocities were omitted with the 
exception of those with respect to w which were retained to 
account for the effect on the horizontal tail of the 
downwash from the wing. If a wing moving through an air mass 


With an angle of attack cdc anda steady forward velocity JU 
O O 
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is suddenly given a vertical velocity change w without 
changing the forward velocity, the angle of attack 


perturbation is aS shown in Figure 1. 





(initial) (Perturbation in w) 


Figure 1. Angle of Attack Due to Vertical Velocity 


Perturbation 


Thus changes in angle of attack, , are proportional to 


the vertical velocity perturbation w. Thus we can write, 





y tw 
O 
tan(m +c SE ———————— BO 
(C8) aes ee re — 
to) O 
and by the small angle assumption 
W 
~ W O 
C= — and © = ( B7 ) 
U O U 
O fe) 
Taking the derivative yields, 
mo Uo ( BS ) 
ro) 
also it can be shown that, 
GZ oe 
= igshiiolee ( BO ) 
Gu w—>0 WwW 
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and substituting equation B 8 yields, 


qZ Pee wl AG 1 

— = Limit —- — =-——Z ( B10 ) 

qw x—-o U & 
re) re) 


also if we redefine the Z dimensional derivatives such that, 


typically, 

; 1 ; 1 P 1 

Z =-—_—_ zZ f Z =— Z y im Z =——Z 

oc U a u U U q ug 
O O 

: ( B11 ) 

Z =e B 

Oe U Pe 
O 


we can now substitute into equation B 3 which yields the 


equation, 


eo-eno + Zu + Z q Sane Zi 6 ( Biz ) 
og u g U Oo e 2 
O 
where Z =1+ Z 
gq gq 


° . ‘ ° ® . 
A Similar expression was obtained for u using the 


Simplifying assumptions and the expressions, 


dX ees | AX. 1 
———- .oeede = ——— = —— X ( 313 ) 
q¥ axa—»~o U oC U oC 

O fe) 


C e ° 
and the expression for u was therefore written, 


uv =X & + Xu+ X'q -g@ cos® + Xy O ( BI4 ) 
co u q fe) oc 
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meres = (xX =— cou } 
q g 0 0 


Finally the expression for q followed directly as, 


6 
——Ne oe ti ut & gq +h 6 (Bio) 
o u g fa 


The model equations for the longitudinal case, including 


an angle of attack due to gusts term, were written in matrix 


form as, 
x = Px + Gu + [w ( B16 ) 
where 
oa (angle of attack perturbation) 
u (longitudinal velocity perturbation) 
state 
x = q (pitch rate perturbation) 
vector 
Q (pitch attitude perturbation) 
oe {angle of attack due to gusts) 
2 
control 
u = Lo. | (horizontal elevator deflection from 
vector e gen @al yi) 
A Te Ze -g/U sine Z 
ao u O 6) oe 
X X a -g cos @ X 
ro u g re) oc 
F = M | M 0 M 
od u g oc 
0 0 1 0 0 
0 —W 
C 


Note: Zz’ in the model is equivalent to the equation of 


motion term (1 + Z) 
g 
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eC 

Ne. 0 

2a) | [w =| 0 
€ 

0 

n 


The measurement equations written in matrix form were, 


y = Hx + Duty ( B17 ) 


eek 
Meee) te =e 
g U 
fe) 
u 
Hx + Du = g v= 
8 





Note: Selected stability derivatives are primed since 
they differ from the notation used in Ref. 6, however they 


are exactly equivalent to the derivatives used in SCIDNT. 
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APPENDIX C 


AIRPLANE AND GROUND DATA SYSTEMS 


mee ALRE LANE DALAT SYSTEM 


1. General 


Se eS =e 


The demonstration F-14A Ship number 8 BU. NO. 157937 
was the test airplane utilized. A detailed description of 
the aircpiane instrumentation system is contained in Refs. 7 
and 9. The airplane was modified to be compatible with the 
Naval Air Test Center Real-Time Telemetry Processing System 
(RTPS). The instrumentation systems and measurenents 
selected for use in the test airplane were primarily 
designed to satisfy the Flight Test requirements for the 
F-14A Performance Evaluation and Demonstration. The test 
airplane was not entirely a production configured F~14A in 
that some miscellaneous cockpit components not reguired for 
flight testing were removed. The production nose radome was 
replaced by a fiberglass and metal radome fitted for an 
instruzentation noseboon assembly. Additionally, some 
production eguipment was reworked in order to accomodate 


various instrumentation subassemblies. 
es Nose Boom Instrumentation 


3 pe ee qe oo ae ce a oe ee eR Ee ee ee ee ee 


The nose boom assembly included Angle of Attack and 
Angle of Side’ Slip Vanes and a Pitot / Static Head. 
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3. Cockpit Instrumentation 


ges eee eee — a 


The Forward Cockpit instrumentation was primarily 
pilot oriented and included; Angle of Attack and Sideslip 
Indicators, Airspeed, Altimeter, Machmeter, Center of 
Gravity Normal Acceleration and an Instrumented Control 
Stick. The Aft Cockpit Instrumentation included left and 
Bg it engine RP i (N1), Main Puel Flow, A/B Fuel 
Flow, Outside Air Temperature and the Instrumentation Master 


Control Console. 


4. Compartment Instrumentation 





~ 


The instrumentation SubcomponentS were mounted 
throughout the airplane in various available locations 
including, the Nose, Forward Compartments, and Equipment 
Compartments. The bulk of the instrumentation equipment was 


located in the forward equipment compartments. 


5. System Description 


= Que a ome 


The instrumentation systen consisted Or the 


following subsystems; 


PCM Signal Conditioning 
PCM 

Timing 

Magnetic Tape Recording 
Telemetry 


Power and Distribution 


O” © -O8O, © OFC 


Ground Support Equipment 


The system had the capacity to record on magnetic 
tape and to telemeter to the ground-data~reduction~ccmplex 
all the PCM data for real-time analysis. The PCM Signal 
Conditioning Subsystem consisted of active pre-sample filter 
Chassis and scme special conditioning required for certain 


measurements. The pre-sample filters had cut-off frequencies 
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of 5, 10, 20, 40, and 80 Hz, The PCM subsystem consisted of 
an analog multiplexer, analog to digital converter, digital 
multiplexer, and the required control logic. The subsysten 
had the capacity of accepting up to 332 analog parameters at 
frequency response from dc to 40 Hz. Discrete data (On/off 
Signals) and a limited amount of digital data and time words 
were time~-division multiplexed into the output data 
strean. The output data was presented as a succession of 
twelve-bit (including parity) parallel data words in NRZM 
format which were recorded on the digital tape. A serial 
output in NRZL format for entry into the "L" Band telemetry 
subsystem was provided. The PCM format is contained in Ref- 
10. An IRIG-£E Time Code Generator was installed to provide a 
common time correlation for magnetic tape and telemetry 
data. The Magnetic Tape Recording Subsystem consisted of a 
PCM Recorder which had the capability to record the PCM data 
and Time Code information. The Digital Recorder contained a 
sixteen track, IRIG compatible head capable of storing one 
hour of digital data on one inch tape at a tape speed of 30 
IPS. The data was recorded in parallel using saturation 
recording technigues. This portion was capable of bit 
densities of 1000 bits per inch. The Telemetry Subsystem 
consisted of a signal combining Und eye an Je Band 
Transmitter and two blade antenna. The operating frequencies 
could be tuned to coinside with those in use at various test 
facilities. The entire Flight Test Instrumentation Systen 
waS powered fron the aircraft primary (generating) 
system. Shouid either the left or right engine or generator 
fail, the instrumentation system would continue until such 
time as both generators failed. The distribution system 
consisted of a 28V dc transformer/frectifier (T/R), and a 
power/fcontrol distribution box, which provided electrical 
system protection and isolation from the aircraft power 
systems through individual circuit breakers. The equiprent 
required to perform checks on the installed instrumentation 


systems was configured in mobile test carts. These carts 
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reguired only 60 cycle power for operation and provided the 
operator with a tool for conducting pre-flight and post 


flight tests on the tape recording systean. 


B. GROUND DATA SYSTEN 


The Naval Air Test Center, Patuxent River, w#aryland, 
employs a xeal-Time Telemetry Processing System (RTPS) as 
universal support for testing on several types of aircraft 
and weapcns systems,Ref. 10. The computer system, designed 
by Xerox, consists of a Sigma 9 general purpose computer as 
the central precessor supported by three Sigma 3 computers 
performing the data input and output processing. RTIPS 
provides the capability to process two complete and 
independent tests (data "streams") simultaneously from 
either tape or telemetry sources. These VStbeans | ean 
Support totally unrelated and asynchronous test flights or 
can support setup for one flight while processing another 
Pine ht . 

Each Data Channel accepts up to 512 channels of PCM data 
and 32 parameters of FM data. This data is processed 
(individually or mixed) by a Sigma 3 and preprocessing 
Programmed Algorithm Unit (PAU). Up to 50,000 values of 
Engineering Units (EU) data per second can be fed to _ the 
Sigma 9 while simyitaneously being stored on digital tapes 
for later playback. Two large screen interactive graphic 
display stations, that provide real-time "hands on" data 
analysis and presentation, are controlled by the Sigma 3 
Display Processor. Batch jobs can be _ processed vhen the 
real-time and non real-time operational modes are not 
utilizing the resources of the system to 100 percent. The 
batch job is defered if additional demands are made by 
either non real-time or real-time tasks during batch job 
processing. 


The specific flight test data utilized by the parameter 
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identificaticn program were, Meee ALPHA FRL, pitch rate, 
pitch attitude, NZ NET, and average horizontal stabilizer 


position. The raw data (digital) was used to generate an EU 
tape for input to the parameter identification program. in 
generating the EU tape corrections were made to the measured 
values as reguired in accordance with Ref 11. For example 
corrections tc angle of attack included nose boom bending 


moments, fuselage bending moments and pitch. 
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